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Abstract 



The a-XY model generalizes, through the introduction of a power-law decaying poten- 
tial, a well studied mean-field hamiltonian model with attractive long-range interac- 
tions. In the a-model, the interaction between classical rotators on a lattice is gauged 
by the exponent a in the couplings decaying as r", where r are distances between 
sites. We review and comment here a few recent results on the static and dynamic 
properties of the a-model. We discuss the appropriate a dependent rescalings that 
map the canonical thermodynamics of the a-model into that of the mean field model. 
We also show that the chaotic properties of the model, studied as a function of a 
display a universal behaviour. 
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1 Introduction 



Equilibrium Classical Statistical Mechanics is one of the grounds of modern 
physics; it originated from the works of Boltzmann and Gibbs at the end of 19th 
century. The computational techniques of Statistical Mechanics, generalized to 
encompass the right enumeration of quantum states, have been successfully 
applied to the study of all states of matter, along all the 20th century. On 
one side, the rigorous analysis of this branch of theoretical physics has been 
impressive and far reaching ||] and the belief that everything in this field has 
been robustly founded seems to be quite widespread. On the other side there has 
been a serendipitous pragmatic use of the Boltzmann-Gibbs approach beyond 
the bounds set by theorems. The use of Equilibrium Statistical Mechanics 
outside the allowed boundaries led to the discovery of anomalies or paradoxes, 
particularly in the field of self-gravitating systems||]Q. In fact all the edifice of 
statistical mechanics rests on few stringent assumptions on the interactions that 
are not fuUfilled by long-range forces. In spite of their fundamental relevance, 
gravitational and Coulombic forces do not fit in with Equilibrium Classical 
Statistical Mechanics (1). 

Short-range interactions guarantee extensivity: that is energy and entropy, 
as functions of intensive internal parameters, grow linearly with the number N 
of microscopic components of the systems. If an extensive system is divided into 
macroscopic parts the total energy and entropy is the sum of the energies and 
entropies of the parts. On the contrary that is not true in systems where long- 
range interactions can reflect themselves in thermodynamic potentials that do 
not scale with the size of the system and cannot be defined in the thermodynamic 
limit (i. e. for N going to infinity). Moreover, long-range forces may well 
induce strong spatial and temporal dynamic correlations that contrast mixing 
and result either in very long relaxation times or in equilibrium distributions 
different from the expected Boltzmann-Gibbs. 

Aiming at a microscopic foundation of non extensive thermodynamics a good 
starting point is the study of the ergodic properties in the thermodynamic limit 
of simple and meaningful Hamiltonian models with long-range interactions Q. 
In the recent past Stefano Ruffo and co-workers have introduced a class of mean- 
field models with infinite range of the interactions . These models are related 
with the statistical physics of self-gravitating systems and that of phase transi- 
tions and became paradigmatic in the study of non-extensivity. In particular, 
the dynamics and thermodynamics of the so called Hamiltonian Mean-Field 
model (HMF) have been extensively studied in the last few years and a compre- 
hensive review can be found in the contribution to these proceedings by Andrea 
Rapisarda and Vito Latora 

The a-XY model we refer to in this paper has been introduced by Celia 
Anteneodo and Constantino Tsallis |^ ; this model cleverly combines the physics 
of the mean-field models quoted above with the statistical physics of lattice 
models of the Ising type with long-range couplings decaying as the inverse power 
a of the distances between sites. Since classical works on the foundations of 
Statistical Mechanics it has been known that these systems are non extensive 
for < a/c? < 1, where d is the dimensionality of the ambient space; a/d is 
then the natural control parameter of non extensivity in this class of models. 

Our work originated from previous works by other researchers active in the 
fields: in the section devoted to the thermodynamics of the a-XY model a 
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canonical solution is presented that has been largely inspired by the analytical 
work of Antoni and RufFo on the HMF model and by a numerical work by 
Tamarit and Antcneodo ||l2| |Q; in the section about dynamics the results on 
the Lyapunov exponents extend the previous work by Anteneodo and Tsallis |^ . 
In the conclusions we try to express our point of view about the connection of 
our results with the particular form of nonextensive thermodynamics proposed 
by Tsalhs 0. 



2 Thermodynamics 

In this section we define the model and we compute its canonical partition func- 
tion. We use the typical methods of Gaussian transformation and saddle point 
integration together with a Fourier diagonalization of the interaction potential. 



2.1 Definition of the model 

We consider the hamiltonian of a system of bidimensional classical rotators (XY 
spins): 

H=i^+y=if:^^^ + -E ''"t'"'^^ - (1) 

The N spins are placed at the sites of a generic d-dimensional lattice, and 
each one is represented by the conjugate canonical pair {Li, 9i), where the i^'s 
are the angular momenta (unit momentum of inertia is assumed) and the ^i's 
€ [0, 27r) are the angles of rotation on a family of parallel planes, each one defined 
at each lattice point. The interaction between rotators i and j decays as the 
inverse of their distance r.y to the power a > 0. Periodic boundary conditions 
are assumed in the form of nearest image convention. 

The hamiltonian is extensive if the thermodynamic limit (TL) — > oo of 
the canonical partition function (In Z)/N exists and is finite. This is assured 
for each a by the presence of the rescaling factor N in front of the double sum 
of the potential energy. A^ is a function of the lattice parameters a,d,N which 
is proportional to the range S of the interaction defined by: 

^«^-E;^- (2) 

The sum is independent of the origin i because of periodic conditions. Then, for 
each a, V is proportional to N . When a > d, which we call here the short-range 
case, S is finite in the TL and things go as if each rotator interacted with 
a finite number of rotators, those within range S. On the contrary when a < d, 
which we consequently call the long-range case, S diverges in the TL and the 
factor 1/A in (|l|) compensates for this. 

The model (H) is intimately connected to its mean-field version for a = 0, 



introduced by Antoni and Ruffo in 1 11 and called HMF model: 
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The model is a classical ferromagnetic XY model and undergoes a ferromag- 
netic transition from a high temperature paramagnetic state crossing the crit- 
ical temperature Tc — 1/2 corresponding to a critical energy Uc — 3/4. The 
computation of the partition function of the full a model, eq. (|l|) is a straight- 
forward generalization of the procedure followed for the HMF model, using the 
standard methods of Gaussian transformation and saddle point technique, the 
whole treatment being successful only for the long-range case a < d. We resume 
here the main passages, whereas a more detailed explanation can be found in 




2.2 Partition function 

The partition function of model (|l]) factorizes in a trivial kinetic contribution 
Zk = (Stt//?)^/^ and a nontrivial potential part 



d^6iexp(-/3y), (4) 



where j3 = 1/{kbT) is the inverse temperature. We rewrite the potential part 
as 

^ = ^. ^-h-}^m. (5) 

where we have introduced an external field h — {hx, hy) and defined a local mag- 
netization vector rrii = {rnix,miy) = {cos 9i, sin 9i). The total magnetization is 
then given by M = {M^, My) = '^iLi ^i- The integral (||) is remanipulated 
introducing matrix conventions. The constraint i ^ j over the double sum is 
removed defining = 1/6, a finite number. Since the numerator l — cos{6i — 6j) 
is zero for i — j the choice of h is free. The removal of the constraint allows us 
to introduce the distance matrix Rij = P/{2Nr'^j). After defining B = /3h and 
C = exp[— /3/(2A^)^-.(l/rf )], the potential part can be written as: 



Zv = C I d^e exp 



(6) 



where fi = x,y. Diagonalizing the distance matrix R — (Rij) we can write the 
first part of the exponent in (H) in a suitable form for application of gaussian 
transformation. The formula requires real and positive eigenvalues: the first 
condition is satisfied because the matrix is symmetric and the second by a proper 
choice of parameter b because, as it can be easily seen, the entire spectrum is 
linearly translated by it. 

Application of the formula leads to known integrals in variables 9i and the 
whole Zv can be rewritten as 



= C'-^ I d''^xd''% (7) 



oo 



,dcti? 

[- ^^.f' + ir i^la^^-f" (|2 *.+B|))] ^ 

where Iq is the zeroth-order modified Bessel function and the 'S'^'s are the gaus- 
sian variables. Use of the saddle point technique is now more delicate than in 
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the mean-field case and requires a careful study of the first and second deriva- 
tive of the function in square brackets at the exponential in eq. (Q) , which we 
call f{w), where w = {^ix, • ■ • , Nx, ^lyj • ■ • i Ny)- The maximum is obtained 
for a vector wo — {"^x, ■ ■ ■ , '^x, '^y, • ■ • , ^y), homogeneous on the lattice sites. 
Defining ^ = {'i>x,'ify), its direction is that of B, and its modulus 5" is given 
by the solution of the self-consistency equation: 

^ = ^{fi[A^ + h]) , (8) 
Jo 



with 



N N 



(9) 



and where Ii is the first order modified Bessel function. We note that when 
h = we have infinitely many degenerate solutions, since only the modulus ^' 
is determined. 

Evaluation of the elements of the hcssian matrix of f{w) at the stationary 
point just proves that it is a maximum. Full application of the saddle-point 
method gives the potential energy contribution to the free energy per particle 
—PF in the TL as f{wo) plus a second-order correction which can be estimated 
as 

^iL- ^EHi-(---)i?k], (10) 

k 

where i?k are the eigenvalues of matrix R. Because of translational invariance 
(in turn due to periodic conditions) they can be computed through Fourier 
transform and thus labelled by the reciprocal lattice vectors k. The N rescaling 
involved in the definition of the matrix R is the key feature to estimate the 
correction. When a < d it is found that most of the eigenvalues go to zero 
in the TL, the quantity being given by N' such that N'/N 0. The same 
estimate then applies to eq. (|l^) and the second-order correction effectively 
vanishes. The procedure cannot be extendend to a > c? because no divergence 
occurs in the spectrum. 

Once obtained the free energy, differentiation with respect to thermodynamic 
variables gives the total magnetization M, which is found to coincide with the 
gaussian variable M = ^P, and the internal energy 

U = — + -(l-M'^)-hM . (11) 



2.3 Universality 

Equations (||) and (|l^) have the exact form of the mean field equations as soon 
as one takes A = 1 or, otherwise stated, as soon as the proportionality in eq. 

is replaced by an equality. As we show below with this position in the 
long-range case the partition function of the model then completely reduces 
to that of the mean-field case and with it the equations of state, namely the 
magnetization curve MvsU, and the caloric curve TvsU where the internal 
energy is U =< H/N >. The thermodynamic of the model thus becomes 
universal in the sense of being independent of the value of a,d for any a < d. 
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As a consequence it shows a ferromagnetic transition at the universal critical 
mean field energy Uc — 3/4. There is a slight difference between an equality in 
eq. (||) and the position A — 1, caused by the term h in eq. but it loses 
significance in the TL because of the divergence of the second term in S. For 
finite A'^ it is a simple matter of computation that the right definition to obtain 
universal state curves is really 

^-S = Y^^. (12) 

Indeed physical equations cannot depend on a pure mathematical parameter, 
such as h. Take for example the high energy form of the caloric equation of 
state. Because the term cos(6'i — 9j) in eq. is null on average (each rotator 
being approximately free), using equipartition of energy we find 2U ~ T + S/ N. 
In order to numerically obtain a universal equation of state we must make this 
expression coincide with the HMF one 2U T+1 and equation ( [T^ ) is required. 

Then we have shown that any model with a < d on any lattice is equivalent 
to HMF if we assume eq. (p^). The theoretical computation supports the simu- 
lative results found in for a onedimensional lattice. As a further check of the 
theoretical result we simulated the system on a simple cubic three-dimensional 
lattice through constant energy molecular dynamics applied to hamiltonian (|l|). 
We used a fourth-order simplectic algorithm with time step 0.02, selected to 
have relative energy fluctuations not exceeding 1/10^. We have chosen a fixed 
N — 343 = 7'^, and have simulated various energy densities H/N and various 
a < 3. In Fig. Q we show that the numerical caloric curves collapse onto the 
universal HMF curve and in Fig. ^ we show the same for the magnetization 
curve. In the same figures we report simulations for one value of a > c? to stress 
the impossibility to extend the previous treatment to the short-range case. The 
curves are qualitatively different: the caloric one has an everywhere continu- 
ous derivative and lacks the finite discontinuity at the critical energy, and the 
magnetization curve has a different critical exponent. As the inset shows if one 
tries to rescale the abscissa for the a > d magnetization curve in order to let it 
coincide with the mean field one, a fundamentally different behaviour appears. 
As known from renormalization group studies in fact the critical exponent for 
the first-neighbour {a = oo) model is ~ 0.36 different from the mean-field 
value 0.5. 

It is also interesting to observe the extrapolation which follows from eq. ([l^ ) . 
If one estimates from this equation a mean value for 1 /r^ one immediately finds 
l/rg ~ N/N. Insertion of this esteem into equation (|^) directly provides the 
fundamental HMF hamiltonian, eq. (||). This remark naturally supports the 
previous analitycal work, however we point out that it cannot fully justify it 
because it leads to a wrong conclusion in the short-range case a > d. 

The physical reason for universality can be given considering a single rotator 
surrounded by all others. It effectively interacts only with the rotators within 
range S (interaction area). For fixed N and varying a this area goes from all 
the lattice for a = to the only first neighbours for a — oo. Increasing N, if 
a < d the range S diverges and the interaction area expands as well, while if 
a > d the range stays finite and the area localized. The idea is shown in Fig. ||. 
When TV ^ oo all interaction areas for a < d simply become all the space, in 
the same way it occurs at finite N only for a = 0. The way the interaction areas 
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diverge obviously depends on a (see note |T5)) but once at the infinity it does 
not matter how you got there. All interaction areas are infinite and equal to 
each other. Simply adapting the interaction field through rcscaling you obtain 
that each rotator feels all the others with the same strength, independently of 
a, d. 

This reasoning cannot be extended to the case a> d. Even for iV going to 
infinity the interaction area remains a well defined localized zone, its features 
depending on a and the conformation of the neighbourhood. 



3 Dynamics 

3.1 Lyapunov exponents in the thermodynamic hmit 

We study the maximal Lyapunov exponent, as defined by the limit p9[ |: 

1 d(i) 



lim X{t) , (13) 



with d{t) = ^/J2ii^^iy + {^LiY being the metric distance calculated from in- 
finitesimal displacements at time t obtained in turn through the double integra- 
tion of both the normal equations of motion 

9, = U (14) 



and the linearized ones 



69, = 5U (16) 

as obtained from hamiltonian (|l|). In the same simulations of the previous 
section we also computed Xmax constructing the curves Xmax{U) for various 
value of a. We report the results in Fig. ^. 

We note a common behaviour for all the curves at a < d. Though not 
perfectly coincident they all show a definite peak around the critical energy Uc — 
0.75, resembling the HMF behaviour as found for example in ||§ [|l) H [|3| || . 
On the contrary the curve for a/d = 5 is different. Besides having a highest 
maximum, it clearly lacks a well defined peak, and it much more resembles the 
first-neighbour behaviour |p5|||2^]. 

We also investigated how Xmax for varying a/d from to oo goes from the 
HMF-like behaviour to the other one typical of the first-neighbour model. 

We simulated the system at fixed energy density Uc — 0.6, varying a/d and 
for N = 343 and TV = 125. The results are shown in Fig. In accordance 
with Fig. ^ the points for a < d are much closer to each other, the vicinity 
increasing for higher N. In the curve for N — 125 a smooth ascent for low a/d 
suddenly increases after a/d = 0.5. The curve for N = 343, though showing a 
similar behaviour, presents the increase for a/d — 0.75. We conjecture that for 
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still higher TV the point of sharp increase moves towards a/d — 1 and then in 
the TL we have 



Kiax{a.,d) = \max{a = <S) ^aKd. (18) 

We checked this conjecture with a more careful study for fixed high energy 
U = 5.0 > 17c- The complete details can be found in j2^. These simulations 
were done using a velocity- Verlet algorithm with a time-step chosen to have 
a relative energy conservation of 10"'* or better. Length of simulations were 
chosen looking at the asymptotic behavior of quantity \{t), eq. (|l3|), where 
d(0) is randomly chosen (see @). We report in Fig. ||, the curves \max{N) for 
various a in dimensions d — Z. 

Because for the HMF model has been theoretically and simulatively 
found that \max{U > Uc, N — > oo) = 0, it is clearly seen that eq. dH) IS 
satisfied. We observe that, for growing N, ii a > d then Xmax is positive and 
constant, whereas if a < c? the maximal Lyapunov exponent tends to zero, the 
value of HMF model. The same behaviour has been found in dimension d = 1 
§ and d = 2 [|7|. 

Then we can believe that in the TL the behaviour of the maximal Lyapunov 
exponent is universal, its dependence on the energy density U being the same 
universal curve which can be found in Fig. 2 of p^ . 

3.2 Finite case 

We have also investigated the scaling law for the reduction with N of the Max- 



imal Lyapunov Exponents. We report here the results shown in |27 . 

For any finite N we have fitted the data with the following functional form: 
Xmax oc N~'^; K being the slope in the log- log plots of Fig. ^. Collecting the 
slopes K as a function of the ratio a/d, for d — 1,2,3, remarkably, through 
the simple scaling a/d, the N dependence of the maximal Lyapunov exponent 
of model (Q) is universal. All curves start at a = from a value close to 
1/3, analytically and numerically found for the HMF model [^|2^; when a/d 
increases from zero to unity all our data collapse into a single curve, and then 
remain zero for a/d greater than unity. And the universal function k, shown 
in Fig. 1^, does not appear to depend on the energy density, provided this is 
greater than the critical one. 



4 Quasi Stationary States 

In literature there is still another interesting feature of the HMF model that 
deserves parallel study and generalization in the a — XY model. Metastable 
states have been observed in microcanonical simulations that give rise to anoma- 
lous thermodynamics and dynamics. In the range of energy density 0.5 < U < 
Uc = 0.75 states with an initial uniform distribution of L^'s (water bag distribu- 
tion) have a slow relaxation and display a negative specific heat [|l) |2^ [|9| |^ . 
The average kinetic energy, anomalously relaxes to a value which correspond 
to a temperature value smaller than the canonical predicted one. The effect is 
purely microcanonical and was also justified theoretically through the micro- 
canonical solution of HMF model . The velocity distribution of these states 
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is not Gaussian g and anomalous diffusion appears in the form of superdiffu- 
sion ||2^ . These states are metastable and called quasi stationary states 

because they have a lifetime which increases with N. The lifetime of these states 
can be estimated either from the relaxation of temperature to its canonical value 
or from the crossover of the mean square displacement from the superdiffusive 
to the normal diffusive behaviour. The two estimates coincide and the lifetime 
is found to grow linearly with N . 

When these states are considered the system then shows different character- 
istics according to the order of limits in the computation of observables. The 
idea is presented in Fig. 5 of . If the thermodynamic limit is carried on be- 
fore temporal average metastable features sussist with their anomalies; on the 
contrary if the iV — > oo limit is taken as the second, normal canonical prediction 
are valid and no superdiffusion is observed. 

We believe these states are present also in the model we presented, eq. (|l|). 
Accordingly we simulated the system at an energy density in the above cited 
range and computed the velocity distribution function at different times, see 
Fig. I^a, b. Though the numerical distribution is constructed only through 
average on particles at fixed time and is not therefore a smooth curve, it does 
not show a relaxation to the Gaussian predicted form (continuous curve) 

The logarithmic scale shows a difference in the tails of the distribution resulting 
in a different microcanonical temperature (variance of the simulative distribu- 
tion) smaller than the canonical predicted one (variance of the theoretical dis- 
tribution). Probably the latter is reached for longer integration times together 
with the gaussian form for the distribution. The affirmation is currently under 
verification. We conjecture that Tsallis' generalized thermostatistics may help 
to describe the features these states show. 



5 Conclusions 

In the first part of this paper we have shown that the a-XY model is canonically 
equivalent to the HMF model. Through a rescaling by N, as defined in (|l^), it 
is possible to map all thermodynamic functions of the a-XY model into that of 
HMF. All a-models have the same universal canonical thermodynamics, that of 
the mean field model for a = 0. In these proceedings a microcanonical solution 
of the HMF model has been given ; a similar solution of the a-model will 
complete the proof of this kind of universality. 

A dynamical equivalence between the HMF model and the a-XY model is 
also suggested by our initial exploration of quasi stationary states in Fig. ||a, 
b. Let us note that what has been said for the thermodynamic limit of these 
states in the HMF model should be extended to the corresponding states of the 
a-model. 

Another kind of universality has been found studying how maximal Lya- 
punov exponents become zero in the thermodynamic limit in the nonextensive 
a < d case. The curve K{a/d) is universal; the way mixing is reduced by the 
range of interaction in these a-models is controlled by the ratio a/d. Firpo has 
given the theory of the scaling of Xmax with for the HMF model |Q . It would 
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be nice to extend this theory to predict the universal form of the K{a/d) curve 
of our Fig. 1^. If the reduction of mixing here reported is connected with Tsalhs' 
nonextensive thermodynamics it should be possible to find out a strict relation- 
ship between Tsallis' entropic index q and a/d; for a/d < 1 and searching at 
equilibrium states in the limit t oo N oo, as discussed by Constantino 
Tsallis in his contribution. 

Friendly and wit remarks by Stefano Ruffo have been very useful while we 
were writing this review paper. 
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Figure 1: The solid line gives the canonical theoretical caloric curve (tempera- 
ture T vs energy density U) for model (|l|) compared with microcanonical sim- 
ulations of the model on a threedimensional simple cubic lattice. Three values 
of a below d (symbols) and one value above (dashed line) are shown. Note the 
qualitative difference between the two cases. 
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Figure 2: The solid line gives the canonical theoretical magnetization curve 

(magnetization M vs energy density U) for the system. The simulative points 
are represented by symbols for the long-range [a < d) cases and by the dashed 
line for the short-range (a > d) case. The difference between the long-range and 
the short-range case is shown in the inset where the impossibility for overlapping 
through abscissa rescaling is evident because of a different critical exponent. 
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Figure 3: Interaction area: it diverges for a < d, is finite for a > d. The long range 
makes all areas infinite and equal to each other. 
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Figure 4: Curves \max{U) ioi d = 3 {N = 343 = 7^) and a/d = 0, 0.5, 0.75, 5. Note 
the similarity between the long-range curves for a < d. 
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Figure 5: \max{a/d) for d = 3, «7 = 0.6, and AT = 343 = 7^, 125 = 5^ Solid curves 
are B — spline interpolations to be used as a guide for the eye. The long-range (a < d) 
points are much closer than the short-range (a > d) ones. The horizontal lines are 
drawn at the simulative values of Xmax found for a/d = oo (first neighbour model). 
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Figure 6: \max vs A*' (log-log plot) for d = 3, = 5 and varying a. Solid lines are 
fits with the functional form (o — -^^ li^Y 
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Figure 7: The mixing weakening exponent k vs a/d for d = 1, 2, 3 (d = 1; from 
d = 2, 3: from [^). It describes the asymptotic A'^ behaviour of the maximal 
Lyapunov exponent \m.ax at fixed energy above the critical one, i.e., Xmax oc 1/N'^. 
The solid line is a guide to the eye consistent with universality. For q = we have 
k(0) = 1/3 (Vd). Figure from 
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Figure 8: Time evolution of probability distribution function (PDF) for the mo- 
mentum p compared with the theoretical canonical one. The graphs refer to a 
three dimensional simple cubic lattice and N = 4096, a = 0.6, U = 0.69. Water 
bag initial (t = 0) distribution is shown in (a) . Microcanonical velocity distribu- 
tion at t = 4000 is shown in (b): the microcanonical temperature defined as the 
variance of this distribution gives the value T = 0.39. continuous curves are the 
canonical distribution computed for T = 0.4757 (see eq. (p^)), corresponding 
to U = 0.69 in the canonical caloric curve, see Fig. Il|. 
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